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SU^f'!ARY 


An aerodynamic analysis system based on potential theory with edge con- 
sideration at subsonic/stpersonic speeds and impact type finite element 
solutions at hypersonic conditions is described. Three-dimensional config- 
urations having multiple non-planar surfaces of arbitrary planform and bodies 
of non -circular contour may be analyzed. Static, rotary, and control longi- 
tudinal and lateral-directional characteristics may be generated. 

The analysis has been implemented cn a time sharing system in ccnjmc- 
tion with an input tablet digitizer and an interactive graphics input/output 
display and editing terminal to maximize its responsiveness to the preliminary 
analysis problem. CDC 175 confutation time of 45 CPU seconds/Mach number at 
subs onic-stper sonic speeds and 1 CPU second/Mach number/attitude at hyper- 
sonic coiditions for a typical simulation indicates that program provides an 
efficient analysis for systematically performing various aerodynamic config- 
uration tradeoff and evaluation studies. 
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INTRODUCTION 


Aerodynamic numerical analysis has developed to a point where evaluation 
of complete aircraft configurations by a single program is possible. Pro- 
grams designed for this purpose in fact currently exist, but are limited in 
scope and aboiaid with subtleties requiring the laser to be highly experienced. 
Many of the difficulties are attributable to the numerical sensitivity of 
the associated solution. In preliminaiy design stages, some degree of appro- 
ximation is acceptable in the interest of modest tum-around time, reduced 
conputational costs, sinplification of input, and stability and generality 
of results. The inportance of short elapsed time stems from the necessity 
to systematically survey a large number of candidate advanced configurations 
or major component geonetric parameters in a timely manner. Modest conputa- 
tional cost allows a greater number of configuratiois and/or conditims to be 
econcmically investigated. 

Qie approach in this spirit is to enplpy panel approximaticns which 
reduce the member of simultaneous equations required to satisfy flow boun- 
dary conditions. Surface chord plane formulations, locally two dimensimal 
crossflow body solutions and non -interfering panel sinplificaticns are exam- 
ples of approximations which can be used for this purpose. 

Finite element analysis when combined with realistic assessment of limi- 
tations and estimated viscous characteristics provides a valuable tool for 
analyzing general aircraft configurations and aerodynamic interactions at 
modest attitiodes for subsonic /supersonic speeds and evalioation of conpressible 
non-linearities at high Mach numbers. 
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LIST OF SYMBOLS 


A Projected oblique cross section area 

Aj_j Influence coefficient. Normalwash at control point i due 

to vortex panel j of unit strength 

Ai Area of quadrilateral panel i 


b 

c 

c 

CAVG 

Cd 

Cd 

CF 

Ci 

C£ 

C£,Cm»Cn 

CL 

Cn 

Cp 

cpnet 


Reference span 
Local chord 
Reference chord 
Average chord 
Section drag coefficient 
Drag coefficient 

Flat plate skin friction coefficient 
Boundary condition for control point i 
Section lift coefficient 

Rolling, pitching and yawing moment coefficients 
Lift coefficient 

Section normal force coefficient 
Pressure coefficient (P-P«>)/q 

Net pressure coefficient CP£*Pu)/<l smd vortex panel strength 


CxjCyjC^ Axial, side, normal force coefficient 
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LIST OF S\TIBOLS (CONTINUED) 


C* 



Uoo 


T» 

T* 


Px»Fy>F2 

gCx) 

h 

K 

i 

£Ci,n) 

L 

L/d 


Force components 

Axisymmetric outer solution to potential equation 
Radius o£ curvature of cross sectional boundary 
Unit vectors in x,y,z direction respectively 

Drag due to lift factor or skin friction thickness correction 
factor 

Equivalent distributed sand gain height or attainable suction fraction 
Effective length 

Length of segment i, i+1 of contour Cn 
Equi\'-alent body length or geometric length 
Body fineness ratio 


M 

n 


p.q,!' 


p 


A A 


p 


Pr 


Mach number 
Moment components 
Unit normal 

Rolling, pitching and yavring velocity about x, y and z 
Nondimens ional angular velocities pb/2U, qc/2U and rb/2U 
Pressure 
Prandtl nxmiber 
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LIST OF SYMBOLS CCOOTIMJBD) 


q 

r 

R 

’'[ 1 
8 . 

s 

S 

Sref 

T 

t/c 

u,v,w 

U 

V 

W 

x,y,z 

x,r,3 

Z 

a 

ai 


Free stream dynamic pressure 1/2dij2 
Recovery factor 

Unit Reynolds number or radius of curvature 
Reynolds number based on [ ] 

Gas constant 
Segment arc length 

Body cross sectional area or surface area 
Referaice area 

Static temperature -°R or tangent of qiiadrilateral panel 
leading edge sweep 

Airfoil thicloiess ratio 

x,y ,2 nondimens ional components of perturbation velocity 
Freestream velocity 
Jet velocity 

Conplex, potential function 
Body axis coordinate system 
C>'lindrical coordinate system 
Complex mmiber y+iz 
Angle of attack 

Local angle of attack at surface control point i 
Angle of sideslip or 
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LIST OF SYMBOLS (CONTINUED) 


y 

V 
5 

n 

ov 

d 

A ■ 
u 

V 
p 

cr 

T 

<P 

$ 

'P 

Q 

Subscript 

c 

CG 

F 

£ 


Vorticity strength per unit length or ratio of specific heats 
Horseshoe vortex strength in Trefftz plane 
Deflection or inpaet angle 

Lateral surface coordinate 
Body slope 

Dihedral angle of quadrilateral panel or boundary layer momentum thiclcness 

Sweep angle 
Absolute viscosity 

Kinematic viscosity, y/p 

Density 

Source density 

Side edge rotation factor 

Perturbation velocity potential 

Total velocity potential 

See figure 3 

Leading edge rotation factor 
camber 

center of gravity 
friction 
lower surface 
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LIST OF SYMBOLS (CONTINUED) 


LE leading edge 

r recovery 

t thickness 

T tip 

TRAN transition point 

u upper surface 

V vortex 

w wave 

“ freestream condition 

Superscripts 

' first derivative or quantity based on effective origin 

" Second derivative 

* Eckert reference tenperature condition 
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SUBSONIC/SUPERSONIC 


The arbitrary configurations which may be treated by the analysis 
are siinulated by a distribution of source and r/ortex singularities. 

Each of these singularities satisfies the linearized snail perturbation 
potential equation of motion 

The singularity strengths are obtained by satisfying the condition 
that the flow is tangent to the local surface : 


All of the resiilting velocities and pressures throughout the flow may be 
obtained when the singularity strengths are known. A configuration 
is composed of bodies, interference shells and aerodynamic surf aces (wings , 
canards, tails etc.). The following types of singularites are \ised to 
represent each. 


wing and vertical tail 



The first step in the solution procedure consists of obtaining the 
strengths of the singularities simulating the fuselage and nacelles, from 
isolated body solution. The- present analysis uses slender-body theory to 





predict the surface and near field properties. The solution is con 5 )osed 
of a compressible axisymmetric con^onent for a body of revolution of the 
same crossectional area and an inconqjressible crossflow conponent, , 
satisfying the local three dimensional boundary' conditions in the (y,z) 
plane. The crossflow is a solution of Laplace's equation 

+ <P =0 

A two-dimensional surface source distribution formulation is used to 
obtain this solution. When the body singularity strengths are determined, 
the perturbation velocities which they induce on the aerodynamic surfaces, 
or other regions of the field, are evaluated. 

The assumptions of thin airfoil theory allow the effects of thickness 
and lift on aerodynamic surfaces to be considered independently. Therefore, 
the effects of the aerodynamic surfaces can be simulated by source and 
vortex singularities accounting for the effects of thickness and lift, 
respectively. The source and vortex distributions used in this program 
are in the form of qioadrilateral panels having a constant source or vortex 
strength. The vortex panels have a system of trailing vorticies extending 
undeflected to downstream infinity. The use of a chordwise linearly 
varying source panel is provided as an option to eliminate singularities 
associated with sonic panel edges at supersonic Mach ntmbers. 

The panels are planar, that is they have no incidence to the free stream 
(although dihedral may be included) , since thin airfoil theory allows the 
transfer of the singularities and boundary conditions to the plane of the 
mean chord. These boundary conditions are satisfied at a single control 
point on each panel. For thickness, the control point is located at the panel 
centroid while the effects of twist, camber, and angle of attack are . 
satisfied at the spanwise centroid of each vortex panel and at 87.5 
percent of its chord. 

A cylindrical, non-circxilar, interference shell, composed entirely of 
vortex panels, is used to account for the interference effects of the 
aerodynamic surfaces on the fuselage and nacelles. The boundary conditions 
on an interference shell are such that the velocity nomal to the shell 
induced by all singularities, except those of the body which it surrounds, 
is zero. The bomdary conditions are satisfied at the usual control points 
for vortex panels. 

The following sections define the details of the solution procedure. 
Incltded are discussions of the isolated body analysis, surface finite 
element analysis considering edge effects, and evaluation of aerodynamic 
characteristics including drag. References are cited for the reader 
interested in further pursuing a particular point. 
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BODY SOLUTION 


According to slender body theory^’2 the flow distirrbance near a 
sufficiently regular three-dimensional body may be represented by a 
perturbation potential of the form 

^ j X) ♦ 5CX) (1) 

is a solution of the 2-D Laplace eqxaation in the y, z cross flow 
plane satisfying the following boundary conditions 

7<f = Jv + kw = 0 

on C(x) (2) 

. A general solution for may be 
written as the real part of a conqilex potential function W(2-) with 2- = y + iz. 

W * JU z £ A_cx) 2* I 

A usefiil alternative representation of and W is obtainable with the aid 
of Green's theorem? 

f ^ W - - 2 (II j <rci) JU.iz-1) Am C3) 

c(x) 

where <rCs) is a "source" density for values of 1 = Vq izcj CYcjZc) being 
coordinates of a point on the contour c(x). 

The function g Cx) is obtained by matching <fi of equation (1) which is 
valid in the neighborhood of the body with an appropriate "outer" solution. 
gCx) is then found to depend explicitly on the Jfech number M and 
longitudinal variation of cross-sectional areas S(x) 


d g 
3-^ 


C(x) and n, are defined in figure 1 


I j s"Cii) Av-(x-f) I- -L j JU.C4-x)^t 

® X. 

- ^'co) ^ .SCO ^(i-x) I 

X 

■^Cx)= I s'cxj - j 3"ct) JU.(x-t) dt J 


M < i 


(4) 


ri > ( 


where 

A • Vl I 
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The body axis perturbation velocities are obtained by differentiation 
of equation (1) 


V 

U5 , 


At st^ersonic speeds , -zone of influence considerations require that u = v = 
w = 0 for X-/&r- < o. 

Solution of the preceding equations is based on an extension of the 
method of reference 3. 

CROSS FLOW COMPONENT 


The reduction of computations to a numerical procedure utilizes the 
integral representation of <$) given in equation ( 3 ) by discretization of 
the cross sectional boundary into a large number of short linear segments 
(figure 2} over each of which the source density is assumed constant at 
a value determined by boundary conditions. 


1 


Computation of «’(i,n) over the segment i, i +1 proceeds by applying ^ 
the boundary condition equation ( 2 ) at each segment of Cn* If =q=Jv+^ 
represents the velocity vector, the corresponding complex velocity in the 
cross flow plane is obtained by differentiation of W in equation ( 3 ) with 
respect to Z: 

V-ItA) « - 2 i (5) 

J A * S 


The contribution by the sources located on segment i, i +1 to the velocity 
at Pj,n is first evaluated. Noting that i, i+1 makes an angle 0(i,n) 
with respect to the horizontal axis, we have 




( 2 . 


and the contribution to the integral in equation ( 5 ) may be written: 




4 


I. 


Z, - S 
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) 


y 


Figure 2. Cross-section Boundary Segmenting Scheme 
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After integration of the last term and summation over all contributing 
segments, the result may be written 






{jU 


R,Ct, j j-w.) 


j (6) 


in which, referring to figure 3 , the quantities R(i,j,n) and S(i,j,n) 
are defined by the relationships 

- 2L. “ ^ . 

4, -A. 




To insure uniqueness of the complex velocity, care must be exercised 
in assigning values to the angles 'p(i,j,n) and 9'(i,j,n). Referring to 
figure 3 , these are measured counter-clockwise from the positive y 

axis so that when facing from P ^ to , a point Pj 

left of i,i+l shall define an angle 'f (i,j,n) = 5(i,n). As Pj ,n 
traverses a path around Pi^n ^ point just to the right of i,i+l, t(i,j,n) 
increases from S(i,n) to ©(i,n) +2ir . The same holds true for f'Ci,j,n) 

as Pj ,n traverses a path around Pi+i,n* consequence of these definitions 
S(i,j,n) becomes -it when approaching i,i+l from the right and rr when 
approaching from the left. This discontinuity reflects that exhibited by 
the stream function upon traversing any closed path which encloses a 
distribution of finite sources. 

From the boundary condition eqmtion C 2 ) , we have 





After substitution of v and w from equation ( 6 ) , this last expression 
becomes 


where 



4,-.^ 




^ o-c. jjt) 



JLw. 


R.(l*I,4 .-w.) 

C i ,4 


+ <■ , j j - 6 t !■ 


( 7 ) 
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p- . 



Figure 3. Details of Variables Pertaining to Segment i,i+l 
of Boundary Cn 
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The sxirface normal perturbation velocity - ^ may be written in terms 

of the body slope Cd‘>'/ax)j^^ , the angles of attack , and sideslip 

and the angular velocities p,q,t" as 







Satisfying equation 7 at each of the points Pj ^ on a given contour 
boundary yields a set of equations for cr (i ,n) . ’ 


AXISYMMETRIC COMPONENT 

Differentiation of g(x) must be carried out with due concern for the 
nature of the improper integrals appearing in eqioation (k) . The result is 


where 


liCx 




3'U,) - ^ I s-vj jU - T. 

L 5'(o)-»- — — -s'ci)- - 5"coJUCi-o\ 

O-x^) J 




>1 


* r* - » 

X^xji j 5 lvv(.X^-t )5 (t)dt = ^ 1 ^ 


’‘-N 




■= 


JUcx^-t> s‘"(tUt = 2. [ ^1*7 1 ^ 


*tC X + 
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To compute the second derivatives of the equivalent body cross sectional 
area required for g'(x3,the first derivatives at % are found by finite 
differences between Xj^ and %+i. Second derivatives S' 'Cm) at xm 
Cm+1 + %i)/^ are then fotind by finite differences between 5' at m and 
%+i. Finally S"Cxm) is determined by linear interpolation of S"(x^) 
between x and S 

PERTURBATION VELOCITIES 

The axial velocity u depends on (a^/»x) and the axisymmetric solution 
g"(x). (3f/3x) is obtained by differentiation of the integral in equation 
( 3 ) to first obtain an exact expression which is then approximated by 
evaluating the result over the segmented boundary . 


The derivation of must take into account the 

path of integration in eqiaation (3) is a function of x. 
figure 1 increments of a dependent variable taken along 
by el( ) and increments taken normal to e are denoted by 
tion of equation (3) then yields 

fact that the 
Referring to 
CCx) are denoted 
SC)* Differentia- 

■|r ' ■ 2 IU| ^ If - j 

<m) y 
'2.-1 

* j <rC5)Ji~(z-S) -^1^ 

j C8) 

From figure i 



C9) 

where h(l) is the radius of curvature of C(x) at ^ . In addition, we have 
from figure i 

Sx ^ 

CIO) 
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To evaluate we note, 

Sx 


Ml 

5x 


S»^o 


<ru j ^*\*4 } ^ ( 4 

_ 


(11) 


Introducing equations (9 ) , (.10), and (11 ) into equation ( 8 ) , 


ax 


^ { i r It] 


Again, assuming that quantities in the brackets of the integrands are 
constant over i,i+l, 

/ \ 5 T J f fhl\ .4. £. lil 1 A 

- <ya,--)(|7). _ SCi,-!,-*-) j 


where 

A<^ tt,j 


•| R.Ct + I.Jj'W) • H t lij-n.) jUv. 


— RC t ,4,-«.) • •n.( Lj-rt.) 

*■ iU t,^) I 
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The radius o£ curvature hCi,n) and the derivatives Sr/Sx, , Sv/Sx ai'e 
approximated at the mid points of the segments i,i+l as follows 

a) h <r/ S X - the derivative at the mid-point x ^ of the interval 

%,xn+l is set equal to the divided difference between <rCi,n) and <rCi>n+l). 
Linear interpolation between these derivatives then yields 5<r/Sx 
at Xn. 


b) S V / S X - referring to figure 4 , the displacement S -y is 

determined by linear interpolation between ^ ^ i and S I i+l,n* 

S "y/Cxn+l - Xn) then represents ^V /Sxat x 'n* ’Linear interpolation between 
the stations x'n then yields Sv/Sx at Xn . 

c) 1/h: - -S at Pi n is detemined by interpolation between values of 
6(i,n) at P i^n- curvature 1/h at P i^n is set equal to the 

divided difference between & at Pi+i n ^ s-'*- i*i ,n • 

The lateral and vertical perturbation velocities, V and , are 
obtained from 


V - 


- 2 


<ni) 

Z-5 




Integration over the boundary with constant segment source density yields: 


V C t cOCj,'*v) 


€ 


<T ( 


16 ( 1 ,-*^) 

e 



R( 1*1 ,j j-yy.) 


L 




} 


Thus 




hd (b 2. ^ I Xw. ^ ^ C t ,«» j'vi) Q 
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SURFACE SOLUTION 


The wing, canard, vertical and horizontal tail are simulated by a 
system o£ swept tapered chord plane source and vortex panels with two edges 
parallel to the free stream. The coordinates of the panel comers are 
specified with respect to an (x,y,z) system having its x axis in the free 
stream direction and its z axis in the lift direction. The panel 
influence equations are written in terms of a coordinate system having a z 
axis normal to the panel and an x axis along one of the two parallel edges. 

A coordinate transformation is necessary to obtain the coordinates in the 
panel reference system. If the plane of the panel is inclinded at an angle 9p 
with respect to the y, z plane, a transformation into the panel coordinate 
system (xp,yp,Zp) is accomplished as follows: 


Xp « X 

*■ 

2 p • ” ©p + a c** ©p 




^ Vp OA(e^-©p> 

T _ Vp >w.(6^-ep> f u)pC*.i(6p.q.) 


A transformation of the (up,Vp,Wp) velocities into the coordinate system of 
the panel on which the control point is located (tic>Vc>Wc) results in the 
axial, binotmal and normal velocities induced on the panel. 

For the image of the influencing panel, the signs of y, 6(- and are 
changed while using the same calculation procedure. 


20 



PANEL SINGULARITY STRENGTHS 


The source singularity strengths may be found directly by equating 
each source panel strength to the slope of the thickness distribution at 
its control point. For panel i 


4 


fAis\ 

Ji 


where refers to the shape of the thickness distribution. The influence 
equations for the source panels can then be used to obtain the velocities 
induced by the source panels anywhere in the flow. 


The determination of the vortex panel singularity strengths are the 
final step in the solution procedure. They are obtained by solving a set 
of simultaneous equations utilizing the vortex panel influence equations to 
relate the singularity strengths to the boundary conditions at the control 
points of the vortex panels. The boundary conditions permit the condition 
of tangential flow to be satisfied. 


Each vortex panel j having singularity strength Cp, induces a set 
of velocities C ^ ) on panel i. Therefore a set of influence 

eqtiations can be written: * 


u. 


C. 

4 



u. 

o. 







4 

£ 

4 




V 



where C^oi, Voi, refer to the velocities induced by all other body 

and source singularities, and written in the coordinate system of the panel 
containing the control point. Since the resultant velocity along the normal 
at a panel control point must be zero. 



and the following system of equations results 



This set of linear equations can be solved for the Cp. and, since it 
assumes symmetrical panel loading, can be used to determine the longitudinal 
characteristics. A similar set of equations exist for the calculation of 
the lateral/directional characteristics. This set assumes an antisyrametrical 
panel loading and has a correspondingly different set of influence 
coefficients Aij . 


BOUNDARY CONDITIONS 

Several types of basic and tmit boundary conditions are considered 
and can be classified as either symmetric or antisymmetric. Linearized 
theory allows the supeiposition of these basic unit solutions. The p, q 
and r rotary derivative boundary conditions are the result of placing the 
configuration at 0, 0 in a flow field rotating at one radian 

per second. 

Symmetric: 

1) basic (<i£e/dL^') - u) - uj, 

B B 


( (Az^cLx ) 


» surface slope due to twist and camber 


a noimalwash induced by slender 
® body thickness and camber 

= normalwash induced by soiirce panels 





2 ) ISiit alpha 


3) Ifeit q rotation 

4) Itoit flap 

Aatisynanetric: 

1} IMit beta 

2) Uhit p rotation 

3) liiit r rotation 

4) IMit flap 


iSo 




« noimalwash induced by slaider body 
at unit alpha 


- y c^&^- U3 


» normalvjash induced by slender bo<fy 
* undergoing unit q rotation 


i8o 


^ « 1. for flap panel 
<r ■ 0. for others 


P‘, 


normalwash Induced by slender bo<^ 
at unit sideslip 




o. 




noimalwash induced by slender body 
undergoing unit p rotation 


b 


Cx-X,^) - U5 

» 


1-3^^ ■ Noimlwash induced by slender bocfy- 
® undergoing unit r rotation 

^ <r 
180 

^ ■ 1. for flap panel 
■ 0. for others 
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CONSTANT SOURCE AND CONSTANT VORTICITY PANEL INFLUENCE EQUATIONS 

The source finite elements have a discontinuity in normal velocity across the 
panel surface while the vortex finite elements have a discontinuity in the 
tangential velocity in a direction normal to the panel leading edge. The magnitude 
of the discontinuity, in each case, is constant over the panel area. In addition 
the vortex panels have a system of trailing vorticies extending undeflected to 
downstream infinity. 

A constant pressure or constant source panel with a quadrilateral shape 
can be constructed by adding or subtracting four semi-infinite triangular 
shaped panels^. These semi -infinite triangles, each determined by a comer 
of the quadrilateral, can be assumed to induce a velocity perturbation every- 
where in the flow. However, each comer represents only an integration limit, 
and all four comers must be included to make any sense. 



If it is kept in mind that four comers must be included, one of these 
triangles having sides determined by y = 0 and x-Ty = 0, induces the 
following perturbation velocities: 
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% 

X ♦ 


/3 


1 - Pt„ 


/S * *') 

t. 


A " 


I fii •» Q 

. i /&*■< o 


constant source panel 


U,(X,^^ 2 ,T) 


q~A [ 

, ft R>x (x-Ta) « fl R + tTx ai . 2f?- 1 

" 4ir 1 


■ x;^-T(^Vt'3l 


' 1 v/t*+a'^ a ♦ CTx-»/a>*a> ] 


+TT [ 

2 " ^/rV‘ R - J 


_ <tA 

+TT 

f 1 (ft R+x T 1 r ♦ Ct» 


1 z R-x z <Vt*V R - CTx*/3.n) 


'tA 

4ir 




constant vorticity panel 


■p « _v 1 I D l^+X ,_i S I l> R ♦ (T»«^i) ,. _ ,4. *• zR 

X = — T” J ^ t"o — — ~ 

8ir [ 2 Ci I^-X ^ ^ a - iT«*ya?ii ^ 

- (z-A.) ^ T 2 Y (.^^*■ 2 ^) + C.J<-T^) * "i* j I 


Vy(^^,2,T) = 

'•’y * 


c.J« 

87 T 


{ 




a iz 


— a-A) -|- 1 


c^A [ 

8 tt [ 


tW 


e R 






zx 


2 




C..A -j_ 9 R*x _ . t i. I !_ 5 -/t\/>^ > 

fiTrl 2^ R-x ■ i C:iV2^) 


25 



The total panel solution is built ip by combining each o£ the four 
comers . 


- j ^ Cx-x^, 2^T*) 

U(X,i,Z) = tx. ( X-X,, , Z,T, ) - T, ) 

u.Cx-Xj, 


etc . 


These restilts hold for both subsonic and supersonic free stream velocities. 
In the latter case, only the real (downstream) contributions are considered. 

The perturbation potential expressions are derived in Appendix A. Sub- 
sequently, verification of the perturbation velocities is presented. 
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LINEARLY VARYING SOURCE PANEL INFLUENCE EQUATIONS 


In supersonic flow constant source panels having a sonic edge have a 
real singularity along an extension of this edge. The singularity occurs 
because : 




e'* Cri/s') -» o 




R, ■t-[TCx-x,)t-A(.a-'3,)3 





« oo 



Control points which are near the extension of this edge will have large 
u and V velocities induced upon than. The singularity can be eliminated by 
using panels which have a source distribution which varies linearly in the 
chordwise direction. The resulting continuous source distribution eliminates 
the singularities. The linearly varying source panel influence equations can 
be found by integrating the constant source panel influence equations with 
respect to x. 


u . = ■ 


lO. 


A 

2TT 

A 

ztr 

in 




_ ^ j_ (S_ R+x . <x-Ta> 1 6 VtVis'r. *CTx+/a\i> _ 

. 71-!^ ■; -TiXZ' - R 

“'‘K i 


These velocity components satisfy the same criteria as the velocity components 
for the constant source panels except that the source strength is proportional to 
x-Ty. The source panel finite elements are constructed with the following 
properties . 
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1. 


All panel leading and trailing edges are at constant (|) , side edges 
are at constant y. 

2. Each source finite element is composed of a pair of chordwise adjacent 
panels . 

3. The source strength varies linearly with chord measured from the 
leading edge of a panel pair, i.e. the maximum value of the source 
strength is proportional to the local chord and attains this maximum 
on the panel edge joining the panel pair. 



The perturbation velocities induced by this panel pair are composed of 
contributions from six comers. 




^ f ' 


If there are N panels in zhe chordwise direction there will be N-1 
singularities or unknora source strengths associated with them, xhe 
linear variation in the source distribution means the value of dc/dx 
must be zero at the leading and trailing edges of each span station. 
This may be an undesired restriction and therefore the use of linearly 
varying source panels is optional. 
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EDGE EFFECTS 


The low pressure created by high velocities around a surface subsonic 
leading edge results in a siaction force. As the edge becomes thinner or the angle 
of attack increases, the flow deviates from potential conditions resulting in a 
progressive loss of theoretical suction and an increase in drag. Generalizing 
a conc^t due to Polhamus^, it is assumed that the leading edge vortex created 
by the detached flow in effect rotates the lost suction force perpendicular to 
the local surface. 

In order to iii5)lement this philosophy, a method of determining the spanwise 
variation of potential suction was developed using linear thin wing theory 
and involves finding the coefficient of the 'Xpc term in the chordwise net 
pressure distribution. The analysis is applicable to multiple surface problems 
of arbitrary planform in the presence of bodies at any Mach number. If the 
chordwise net pressure distribution on a thin wing at any given span station 
is expanded in a series 
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it is shown in appendix B that the leading edge nondimens ionalized suction force 
per unit length is 
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where T = 

Z « 

and c is the local chord 


(13) 


Only the first term in equation 12 contributes to the thrust, since it is 
the only contribution which is infinite at the leading edge. If the chordwise 
pressures are known at M points along the chord, the coefficients Aji are 
obtained by fitting a least square error curve described by N terms of the 
series, through the points, where N < M. The pressure distribution is 
obtained using constant pressure panel analysis. 
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The method used to confute the suction force at surface tips is similar 
to that for the leading edge. By using the .irrotational property of the 
flow, it is shown in appendix B that the tip suction force is: 
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The sectional leading edge suction attained in the real flow, Cg (y) , is 
estimated by^ 
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The chord of the normal section, , is defined so as to place the maximum 
thickness , t at the mid chord as indicated in the following sketch. The 
associated leading edge radius is designated by 


Leading 



Potential tip suction is asstmed to be fully rotated as a result of vortex 
formation in the present analysis. 
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AERDDWMIC CHARACTERISTICS 


Longitiadinal and lateral -directional forces and moments due to thick- 
ness, twist and camber, pitch, sideslip, and the dimensionless rotary 
velocities and ^ are obtained from surface pressure integrations of the 

various configuration conponents . 

BODIES 

The pressure coefficient, to an approximation consistant with slender 
body theoiy, is 
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The forces and moments are obtained frcmi the surface integrations 
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In terms of these expressions 
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Crosscoipling between the pitch, sidesHp, and rotarv’ motions through 
the' product and quadratic terms in equation (20) is neglected. 


PLANAR COMPONENTS 


Surface pressure distributions are calculated for planar components 
using the first-order linearized form 








4 


] 


The +/- signs refer to the upper and lower surfaces respectively. The term 
consists of the velocities induced by the isolated bodies and 
other vortex and source panels. These velocities are obtained by multiply- 
ing the ^ influence matrices by the appropriate panel strengths, fhe 
term accounts for the perturbation velocity induced by the 
local distribution of vorticity and changes sign from t^jper to lower 
surface. The total and values are the result of taking linear 

combinations of all the basic and unit solutions. 

The net pressures for each of the basic and unit solutions are integrated 
numerically to give the section forces and moments, component forces and 
moments, and configuration forces and moments . 

Since the vortex panels have a constant pressure distribution, a block 
integration scheme is employed. With the exception of drag, these basic 
and unit force and moment coefficients are combined in a linear manner to 
produce the aerodynamic characteristics for any desired flight condition. 

Since drag varies in a parabolic manner, it must be considered on a point 
by point basis as defined in a later section. 

The longitudinal normal force distribution on the bodies is calculated 
for each solution. The load distribution on the interference shell portion 
of the body is given by integrating over all vortex panels at a given 
longitudinal station. 

normal force 
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where N is the number of panels around the shell, L is the length of the 
body,^X is the length of the interference shell segment, A, is the panel 
area, and Ci = 2 for a centerline body or Cy = 1 for an off centerline body. 

This carryover load distribution is added to the previously calculated 
isolated body longitudinal load distribution. 

The section characteristics of planar components are determined by a 
chordwise sunmation of panel data at each span station and are given by 
the following equations: 


local normal force 
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where N^is the number of chordwise panels, and As is the width of the span 
station and is given by 
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Caimonent forces and moments including edge vortex effects are given 
by the following equations: 
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where N is the number of vortex panels on half of a symmetrical conponent (or total 
for an asymmetrical conponent) and FI, F2 are given by 


symmetric loading 


antisymmetric loading 


Fi = 1 asymmetric geometry 
= 2 symmetric geometry 
F2 = 1 asymmetric geometry 
= 0 symmetric geometry 

Fx = 1 asymmetric geometry 
= 0 symmetric geometry 
F 2 = 1 asymmetric geometry 
* 2 symmetric geometry 


For the leading and side edge vortex terras, Ng is the total mmiber of spanwise 
panels for both component halves, Nct is the number of tip chordt^ ise pan els, XycP 
is the axial location of tip vortex center of pressure, Ag' = AS / 1 + and the 
rotation factors and •f’ are derived in appendix B and defined below. 

Leading edge vortex rotation: 
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where 6 is the slope angle of the leading edge camber line and the sign of coefficient 
A o (from equation 13) is used to determine the direction of vortex rotation. 
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Side edge vortex rotation; 
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where 5 is the slope angle of the tip camber line, ± is plus for the left side and 
negative for the right side of the configuration and the sign of coefficient Cn^ 
(from equation 14) is used to determine the direction of vortex rotation. 


The X coordinate of the center of pressure is given by 
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For interference shell components , the total forces and moments of the 
corresponding isolated body are added to those of the shell. 

The forces and moments for the conplete configuration are obtained by summing 
those of the indiTfJ.dual conponents . 
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DRAG ANALYSIS 


Estimation of configuration aerodynamic efficiency requires the 
calculation of drag. The analysis separates the computation into skin- 
friction and pressure drag coirponents that are assumed to be independent 
of each other. The following form is considered and produces non- 
parabolic polars as a result of the incorporation of attainable suction 
considerations . 
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The specific techniques used for the various drag evaluations are discussed 
below. 


SKIN FRICTION 


Several well established semiaitpirical techniques for the evaluation 
of adiabatic laminar and turbulent flat plate skin friction at incompressible 
and compressible speeds are used to estimate the viscous drag of advanced 
aircraft using a conponent buildup approach.. A specified transition point 
calculation option is provided in conjunction with a matching of the 
momentxM thickness to link the two boundary layer states. For the 
turbulent condition, the increase in drag due to distributed surface 
roughness is treated using unfomily distributed sand grain results. 

Component thickness effects are approximated using e::!q3erimental data 
correlations for two-dimensional airfoil sections and bodies of revolution. 

Considerations such as separation, component interference, and discrete 
protuberances (e.g. antennas, drains, aft facing steps, etc.) must be 
accounted for separately if present. 
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In the following, a discussion is presented for a single corrponent 
evaluation in order to simplify writing of the equations and eliminate 
multiple subscripting. The total result is obtained by a surface area 
weighted summation of the various conponent analyses as described on 
page 65. 

Laminar/Transition 


A specified transition option is provided in the program. The principal 
function of the calctilation is to provide the conditions required to 
initialize the ttirbulent solution. In particular , the transition point 
length and momentum thickness Re>Tiolds numbers are required. 
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This solution is based on the laminar Blasius result C8» chapter \n;i} 
in conjunction with Eckert's compressibility transformation^. This 
option permits an assessment of the reduction in skin-friction drag if 
laminar flow can be maintained for the specified extent. It does not 
establish the liklihood that such a condition will be realized in practice 
or to what extent. 
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Turbulent 


Smooth and distributed rough surface options have been provided in the 
analysis. In either case, the solution is initialized by matching the 
momentum thickness at the transition point produced by the laminar/ transit ion 
solution. That is, an effective origin (commonly referred to as a virtiial 
origin) is established for the turbulent analysis. 

For the hydraulically smooth case 
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where 
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The compressible turbulent flat plate method used here is that proposed 
by Van Driest based on the Von Karman mixing length hypothesis in 
conjunction with the Squire-Young formulation for profile drag (8 , chapter 
XXIV) as applied to a flat plate. 
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For the distributed rough case 
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The turbulent flat plate method ■used here is that of Schlicting 
(8, chapter XXI) which is based on a transposition of Nikuradse's densely 
packed sand grain roughened pipe data. The effect of congress ibility is 
due to the reduction in density at the wall as proposed by Goddardl^ 

The selection of the equivalent sand grain roughness for a given manufac- 
turing surface finish is made with the aid of Table II which was taken from 
Clutter ! ^ 
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TABLE II 


Type of Surface 

Equivalent Sand Roughness 
K_s Cinches) 

Aerodynamically smooth 

0 

Polished metal or wood 

0.02 - 0.08 X 10’^ 

Natural sheet metal 

0.16 X 10“^ 

Smooth matte paint, carefully applied 

0.25 X 10‘2 

Standard camouflage paint, average application 

0.40 X 10"^ 

Camouflage paint, mass -production spray 

1.20 X 10'^ 

Dip-galvanized metal surface 

6 X 10’^ 

Natural surface of cast iron 

10 X 10-3 


Thickness Corrections 


The foregoing evaluations produce as estimate of the shearing forces 
on a flat plate (at zero angle of attack) for a variety of conditions. As 
an actual aircraft has a non vanishing thickness, an estimate of pressure 
gradient effects on skin friction and boundary layer displacement pressure 
drag losses is required. A common procedure for accomplishing this and the 
one which will be used here is based on non- lifting experimental correlations 
for symmetric two-dimensional airfoils and axisyrametric bodies. The 
following relations derived by Homer ( 13 , chapter VI) are used* respectively. 
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Homer reconiinends K, = 2 for airfoils with naxiiTiuin thickness at ZQ-s 
chord and Ki= 1.2 for MA.CA. 64 and 65 series airfoils. In this regard, the 
best information available to an analyst for his particular contour should 
be used. This is especially true for modern high performance shapes such 
as the supercritical airfoil, etc. 

Total Viscous Drag 


The aircraft total viscous- drag coefficient is estimated by a sum of 
the preceding analysis over all components (i.e. wing, fuselage, vertical 
tail , etc . ) . That is 


Co = 


I 



The component length used in the calculation of the skin friction 
coefficient is the mean chord for planar component segments and the physical 
length for bodies and nacelles. 


BASE DRAG 

Blunt base increments are estimated at subsonic and supersonic speeds by 
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The expressions for the base pressure coefficient are derived from correlation 
of flight test results for the X-15, various lifting bodies, and the space 
shuttle. Power effects are treated as reductions in base area in the present 
analysis . 
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POTENTIAL DRAG 


One hundred percent suction drag due to lift and siroersonic wave drag 
due to thickness can be evaluated by integration of the raomentura flux 
through a large circular cylinder centered on the x axis and whose radius 
approaches infinity. 
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The resulting ei-cpression for the total pressure drag is as follows: 
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Tne first tena 'represents tlie wave drag due to Tiomentura losses thru 
the side of the q/linder caused by standing pressure '^aves. The second 
term represents the vortex drag which arises from the kinetic energy/ left 
behind in th.e Trefftz plane by the system of trailing vorticies. 
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Vortex Drag 


The vortex drag may be con5)uted when the distribution of trailing 
vorticity in the Trefftz plane is known. The assim^jtions of linearized thin 
wing theory result in a vortex sheet which extends directly downstream of 
all lifting surfaces. By changing a surface integral for kinetic energy to 
a line integral over the vortex sheet in the Trefftz plane the following 
expressions for lift and drag result . 
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where 


C vortex sheet branches 

weighted section normal force coefficient a. 

^eo asyrptotic normal velocity on the vortex sheet 
vortex sheet branch coordinate 

® inclination of vortex sheet with respect to y axis 
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The program conqjutes the normal velocity on the vortex sheet, , 

by assuming the vortex sheet is congjosed o£ finite trailing horseshoe 
vorticies whose strength is proportional to the local section CnCs) . The 
noimal velocity is confuted at a control point located midway between the 
trailing vortex segments. 



Each vortex j induces a contribution to the noimal velocity at 
section i. 
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The integral for wave drag 


A. 

may be simplified by allowing the cylindrical surface of integration to 
recede infinitely far from the disturbance. Under these conditions, the 
spatial singularity simulations can be reduced to a series of one- 
dimensional distributions. The basis for this reduction is the finding by 
Hayes (1^) that the potential and the gradients of interest induced by a 
singularity along an arbitrary trace on a distant control surface, say PP’ 
of figure 6 (or alternately described by the cylindrical angle , is 
invariant to a finite translation along the surface of a hyperboloid 
emanating from the trace and passing through the singularity. As the apex 
of the hyperboloid is a great distance away, the aforementioned movement is 
along a surface which is essentially plane; it will be henceforth referred 
to as an "oblique plane". Since a singularity is a solution of a linear 
differential equation, all singular solutions which lie on the surface of the 
same hyperboloid (oblique plane) may thus be grouped to fom a single equi- 
valent point singularity whose strength is equal to the algebraic sum of the 
individual strengths and which induces the same potential (momentum) along 
the trace as the group of individual singularities. 

This finding provides the basic technique for reducing a general spatial 
distribution of singularities to a series of equivalent lineal distributions. 
This is accomplished by surveying the three-dimensional distribution longi- 
tudinally at a series of fixed cylindrical angles, © . At each angle, the 
survey produces an equivalent lineal distribution by systematically cutting 
the spatial distribution at a series of longitudinal stations along its 
length. At each cut, the group of intercepted singularities is collapsed 
along the "oblique plane" to form one of the equivalent point singularities 
con^rising the lineal distribution . 
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Tlie £ar-£ield e^q^ression. for the wave drag of a general system 
of lift, and side force elements is 
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is the equivalent lineal singularity strength at the 
cylindrical angle & 

•fCt,©) =! equivalent source strength per unit length 


■T d 3 («,©) = equivalent lifting element strength per unit length 
~ U 3 C<,a) = equivalent side force strength per unit length 
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These strengths are deduced from the three dimensional singularity 
distributions by application of the superposition principle along 
equipotential surfaces. For a distant observer such surfaces are planar in 
the vicinity of the singularity configuration. The individual singularity 
strengths are related to the object under consideration by the requiranent 
of flow tangency at the solid boundary. Lomax [1 5) derived the following 
approximate expressions between the equivalent singxolarity strengths and a 
slender lifting object. 
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where (see figure 7 ) 

A (e,©) is the Y-Z projection of the obliquely cut 

crossectional area 

c is the contour around the surface in the oblique cut 
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Figure 7 Areas and Forces Pertinent to the Evaluation of 
Wave Drag from the Far Field Point of View 


Utilizing the singxilarity strength expressions derived by Lomax, the 
following expression for wave resistance based on the far field theory of 
Hayes is obtained 
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In order to facilitate subsequent discussion, the above result is 
manipulated into the following form 
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A requirement for this transformation is that 
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In accordance with equation 22, the wave drag of a configuration is 
the average of the wave drag of a series of equivalent bodies of revoluation. 
The drag of each of these bodies is calculated from a knowledge of its 
longitudinal distribution of normal cross-sectional area. For each 
equivalent body, these areas are defined to be the frontal projection of the 
areas and the accumulation of pressure force in the theta direction 
intercepted on the original configuration by a system of parallel oblique 
planes each inclined at the given ^iach angle. The common trace angle (<b) 
of the system identifies the equivalent body under consideration. 

Nacelles are assumed to swallow air supersonically. That is, the duct 
is operating at a mass flow ratio of unity. Consistent with this asstm^Jtion, 
the equivalent body cross sectional area distribution is increased by the 
oblique projected duct capture area at all stations ahead of the duct 
which are intercepted by an oblique plane. 

Blunt base components are extended Cmaintaining constant cross 
sectional area) sufficiently far downstream to prevent flow closure arotind 
the base. 

In addition to a geometric description, a definition of the pressure 
distribution acting on the configuration is required. The vortex panel 
analysis is used for this purpose. The thickness pressures for planar 
components have tacitly been neglected under the assumption that the surfaces 
are sufficiently thin that the net pressure coefficient is representative of 
pressure acting on the oblique section. 

Estimation of the wave drag based on equation 23 depends on solution 
of integrals of the type 
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of a numerically given function G(X). Evaluation of such forms has been 
studied by Eminton ^6,17 for functions having G ' CX) continuous on the 
interval (0,1) and G' (o) = G’(l) = 0. In such situations ,G ’ (X) can be ■ 
expanded in a Fourier sine series. It can then be shown that 
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where 
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Eminton then solved for the value of the Fourier coefficients which restolt 
in I being a minijnum, subject to the condition that the resulting series 
for G(X) be exact for an arbitrarily specified set of points o,l, ^i, i = Ij-n.. 

This approach produces the following result 


r . f C;C. f,. 




where 


c. •= 

4 

<; - [aco 




ZC 1-26;) V 6. (l-€;) 

€. * 

— i— 1 

6 t 6 'v^. 

4 



= 

{ " 



+ Z. €.€. ) y 6i«4 0-€;)0-£i) 


The solution of equation 23 for wave drag is accomplished by use of the 
following identities. 
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ZERO SUCTION DRAG 


The zero suction drag due to lift is calculated by numerically 
integrating the net pressure distribution times the projected area in the 
streamwise direction over each of the planer surfaces. The following 
block integration scheme is used to sum over all quadrilateral panels. 
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and 




is due to t'^^^LSt and camber, S is the control surface deflection, 
and dfs 1 for control surface panels and ’5"= 0 for non control surface 
panels. = 2 for symmetric geometries and Fi = 1 for asymmetric 
geometries . 


Edge forces are not considered in this evaluation. 
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DRAG DUE TO LIFT 


The drag due to lift for the total configuration is based on linearized 
potential (100 percent leading edge suction) calculations plus corrections to 
account for suction losses and associated edge vortex forces 
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and the leading edge and side edge rotation factors and Tp, are derived 
in Appendix B and defined below. 
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where 5 is the slope angle of the camber line perpendicxiLar to the leading edge 
and the sign of coefficient Aq (from equation 13) is used to determined the 
direction of vortex rotation. 
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where 6 is the chordwise slope angle of the tip camber line, plus refers to 
the left side and negative to the right side of the configuration and the sign 
of coefficient Gno* (from equation 14) is used to determine the direction 
of vortex rotation. 

As estimate of the average level of leading edge suction for the complete 
configuration is based on the following eqimtion: 

CDt ■ CDl 

SUCTION = ? where Cdt = Cpr for Kc = 0 

" %00 ^ 
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HYPERSONIC 


High Mach mmber analysis is based on non- interfering constant pressure 
finite element analysis^^. 

An arbitrary configuration is approximated by a system of plane 
quadrilateral panels as indicated in the following sketch. 




The pressure acting on each panel of a vehicle component is evaluated by a 
specified congress ion- expans ion method selected from the following Table. 
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Impact Flow 


Shadow Flow 


L Modified Newtoniaxi 

2. Modified Newtonian+Prandtl- Meyer 

3. Tangent wedge 

4. Tangent -wedge empirical 

5. Tangent-cone empirical 

6. OSU blunt body empirical 

7. Van Dyke Unified 

8. Blunt -body shear force 

9. Shock-expansion 

10. Free molecular flow 

11. Input pressure coefficient 

12. Hankey flat-surface empirical 

13. Delta wing empirical 

14. Dahlem-Buck empirical 

15. Blast wave 

16. Modified tangent- cone 


1. Newtonian (Cp = 0) 

2. Modified Newtonian+Prandtl-Meyer 

3. Prandtl-Meyer from free-stream 

4. OSU bl\int body empirical 

5. Van Dyke Unified 

6. High Mach base pressure 

7. Shock-expansion 

8. Input pressure coefficient 

9. Free molecular flow 


A discussion o£ the various methods is presented in appendix C . Specific analysis 
recommendations are provided by the program on a conponent by conponent basis. 

In each method, the only geometric parameter required for determining panel 
pressure is the inpact angle, 6, that the quadrilateral makes with the free-stream 
flow or the change in angle of a panel from a previous point where 
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Voo = Voo cos a cos B i - v„ sin B j + V* sin a cos B k 


Q=pi-qj -rk 


r = (x-xcg) i + Cy-Ycg] j + (z-Zcg) ^ 


Panel switching between it^jact or shadow conditions is based on 6 > 0 in the former 
case and 6 < 0 in the latter. 

AERODYNAMIC CHARACTERISTICS 


The presstire on each panel is calculated independent of all other panels (except 
the shock-expansion method). If the vehicle is rotating, the local pressure coefficient 
must be corrected back to free -stream conditions. That is 
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Vehicle component forces are obtained by sxmiming panel forces 
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■ Sref[ 2® cy-ycgi "z 


ACm= i r 

cSREpL 


aQi 


where 

A = panel area 


x,y,z = panel centroid 


Configuration buildup and total vehicle coefficients are obtained by appropriate 
summation of conponent contributions. 

The conversion of the body axis force characteristics to lift and drag coefficients 
is based on the standard trigonometric relations. 
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The vehicle static stability derivatives in angle o£ attack and sideslip are 
calcialated by the method of small perturbations. Since the basic force and 
moment characteristics are non-linear^ these parameter vary with attitude angle 
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The damping derivatives due to vehicle rotation rate are given in a 
similar manner 
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The control surface derivatives are also calculated by the method of small 
perturbations. 
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CCNCLUSICNS 


An aerodynamic configuration evaluation program has been develqjed and 
implemented on a time sharing system vdth an interactive graphics terminal 
to maximize respcmsiveness to the preliminary analysis problem. 

The solution is based on potential theory with edge considerations at 
subsonic/sL 5 )ersonic speeds and impact type finite element analysis at hyper- 
sonic conditions. Three-dimensional configurations having multiple non -planar 
surfaces of arbitrary planform and bodies of non-circular contour may be 
analyzed. Static, rotary, and control longitudinal and lateral-directional 
characteristics may be generated. 

CDC 175 confutation time of 45 CPU seconds/Mach number at subsonic-stfer- 
sonic speeds and 1 CPU second/^fech number/ attitude at hypersonic conditions 
for a typical simulation indicates the program provides an efficient analysis 
for systematically performing various aerodynamic configuration tradeoff and 
evaluation studies. 
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APPENDIX A 


SUBSONIC/SUPERSONIC FINITE ELEMENT DERIVATIONS 


VELOCITY PERTURBATION POTENTIAI. 


Tne expressions for the velocity potential induced by the source and vortex 
finite elements can be derived from the velocity potential for a point source! 

Source at ^ 




A.q- 
+1T R. 




Therefore the velocity potential for an area in the Zq = 0 plane hnving 
constant source density is: 
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A doublet at Xq, iJo> ^o* 


is the derivative of a point source: 
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Integrating from to infinity yields the potential for a line 

doublet or elementry horseshoe vortex. 
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And an area of constant vortex strength is obtained by integrating this 
expression over the panel area: 


2) 




— — I (2-a) ^ I 



The solution to these integrals is performed in the following sections. Ml 
integrals may be checked using tables 1 and 2 at the end of the Appendix. 

The velocity expressions may be obtained by differentiating the velocity 
potentials using table 1. 
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SOURCE PANELS 


First the integration is perfomed over the panels in the Xo 
direction. 



To integrate with respect to y, a change of variables is introduced: 
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which is independent of /a 
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Hierefore the integral may be performed with i =» constant 
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which may be checked using table 2. 


Each of the four integration liniits corresponds to a comer of the 
quadrilateral. Placing the origin of the Xo, y© coordinate system at one 
comer, the contribution to § becomes: 
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and combining each of the four comers: 
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VORTEX PANELS 


Analogous to the source panels the integration is first performed in 
the Xo direction. 
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changing variables and integrating with respect to 7 
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therefore for one comer or integration limit 
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VERIFICATION OF THE PERTURBATION VELOCITY EXPRESSIONS 


To establish that these are the correct pertxirbation velocities the 
fallowing criteria must he met: 


1. laplace's equation must he satisfied 




■ ** * a:\ 

or the equivalent. 


2 « 


U. « V 
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2, The correct discontinuity or jTnnp in the perturbation velocity must 
occtir at the surface of the quadrilateral panel area* For the source 
■ panel the jimg) occurs in the normal or w velocity and on the vortex 
panel there must he a jump of constant magnitude in the u perturbation 
velocity over the panel area. The perturbation velocities shotild be 
continuous elsewhere, except on the trailing vortex sheet of the 

vortex panel. 

% 


3 . The perturbation velocities must go to zero as upstream infinity is 
approached . 


i*-. For the vortex panel the trailing vorticity must extend straight bach 
to downstream infinity. This means that any discontinuity in the v 
velocity must be zero outside the spanwise bouiaiaries of the panel ard 
must be zero upstream of the panel. 
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The first criteria can be established by using the derivatives given in 
table 1 . 


The second criteria can be established by noting that all terms except 






and 




are continuous at Z = 0. Consider these terms keeping in mind th^t the 
contributions from all four comers must included. 


If we let 
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then the contributions from both comers on the leading edge can be combined 
as follows. 
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I£ we define 
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Therefore when a similar procedure is carried out for the trailing edge of 
a source panel we obtain the following jump in the w perturbation velocity. 
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For the vortex panel (subsonic) vre have an additional term. Considering 
both additional terms from the leading edge comers : 
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The contribution from each panel comer is: 
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Therefore summing all four panel comers 
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To verify the third criteria we must show that all of the functions 
approach zero when all four comers are considered as x - ea 


= - 2- Oi-y(.y * -j /2> 


Therefore considering both corners on the leading edge of the panel 


■* Ixl 

X, 


ii:. 



R, » Cx*x.) 
R, Cx-x, ) 


I 

t 



R,,* Cx-Xt'^ 
Rj- Cx-x^) 





R, f (,x-x.> 
2 * 


O.V 



R ♦ CTx^/aa^ 
R. - CTx +/6*:a^ 


z 



+ T 


and therefore this limit is also zero when both comers of the leading or 
trailing edges are considered. Since all terms are accoxjnted for, the 
perturbation velocities are zero far upstream. 
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there is an apparent singularity along the line 
Cx-T^i^ » O ^ 2. » O 


However this singularity nay be removed by combining the contributions from 
both comers of the leading or trailing edges of the panel. Along either of 
these edges the values of 


Ot- x.^ - TCtk--i. ) and a 
are the same for each of the panel corners. 



It can be seen from the above diagram that CTx will have the same sign 

on a point (x,y,o) which lies outside the spanwise boundaries of the 
quadrilateral. Therefore outside the spanwise boundaries the term 

can be canceled by combining both corners, and the resulting term 
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will not be singular if the correct + or - sign is chosen. >7ithin the spanwise 
boundary an actvial singularity occiirs on the panel edge. 
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a + X 


also has a possible singularity. This term 


The term 
can be 'vnritten 


i log 
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For the source panel the singularity may be removed for points along 
which are outside of the panel boundaires. 


If ( x-x,) and ( X -x^) have the same sign the combination of the two terms 
gives 
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where the correct sign is chosen to remove the singularity. On the panel 
edge the singularity is real and cannot be removed. 



removable singularity 


real singularity 


For a vortex panel the terms (subsonic) 
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Both have real singularities for x >o (downstream) and removable singularities 
for x<o (upstream). The real singularities occizr on the panel edges and on 
the edge of the trailing vortex sheet. 


if X < o 
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SUPERSONIC VELOCITIES - SPECIAL CONSIDERATIONS 


The velocity perturbation influence equations for supersonic flows are 
treated by talcing only the real parts of the expressions. This means that 
R. -s y is set eqiial to zero for points which lie outside the 
downstream f^ch cone from any given comer. Therefore, R and log ■ are 
zero for points which lie outside the downstream tfech cone. For >e>, 

there are no problems using this method. 


If 
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the real part of 
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therefore, combining two comers 
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If a = 0 and either R, or R^ is zero and we allow the other to approach zero, 
the value of becomes 
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Therefore if R, and R^ are zero but we are inside the envelope of mach cones 
from the leading edge (see figure 5)^ the value of F2 is set equal to 
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Figure 5 Siroersonic Leading Edge Mach Cone Envelope 
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APPENDIX B 


SURFACE EDGE FORCES 
LEADING EDGE POTENTIAL SUCTION 

In the limit as the wing thickness goes to zero, the increasingly redijced 
pressure acting over a decreasing area results in a liraitii^ suction force at the 
leadii^ edge. If we consider the leading edge region, the force on the 
airfoil may be obtained by integrating over a control surface in the flow, 
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where S is a control surface into vdiich the leading edge penetrates and F is the 
force on the area enclosed by S. In tvro dimensions the surface integral becomes 
a line integral and since for incon^ressible, irrotational flow 
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viiere C is ti» contour arwaid the leading edge of the airfoil and is the force 
per unit of leading edge length. 
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As the wing thidmess approaches zero, the wing becomes a line segment and 
the flow in the leading edge region is identical to the flow around a 180 degree 
comer. Inconpressibly , it is described by 
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where (r,9) is a coordinate system centered at the leading edge, 
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To relate -Fx, the leading edge suction force, to the pressure distributim near the 
leading edge, the ACp across the line segment must be evaluated. 
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and if c is the chord length 
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These expressions relate the leading edge thrust coefficient to the net distribution, 
ACp, at the leading edge. 
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In general we can write 
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Ao is the coefficient of the y term and therefore determines the leading 
edge suction force since only the term which is infinite at the leading edge 
contributes to the suction. 

c » “X ^ A * o 

For linearized con^nressible flow the following Mach mmiber correction must be applied 
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Let the subscript or supercript o denote tfie variables normal to the leading 
edge. Then 


^ Cp ■* Pi AiVn. ^vv.^ 


and 
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the ratio of thrust per unit length is identical in either system 
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ACp and C-t in the freestream coordinate system are based on freestream dynamic 
pressure q<». Thus 
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therefore combining terms 
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when ACp is given by 


ACp ^ f^o ^ if 2. 


SIDE EDGE POTENTIAL SUCTION 

The method used to con5)ute the suction force at surface tips is similar 
to that used for the leading edge. Since the flow is irrotational 
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introduce a change o£ coordinates 


let S be the fraction of chord 

T be the slope of a constant 5 line t(s_, 7 ) 



then integrating 
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Near the tip, we assume a net pressure coefficient of the form 


aC^C5,7)= — - ^ L , 


where 
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Differentiating 
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Then as n -> rtmax* equation 2 gives (keeping only the largest term) 
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Using the expression derived for flow around a comer (equation 1) in conjunction 
with this relation, the suction force at the tip is given by 
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EDGE FORCE AND M»ENT INCRBENTS 


To account for edge vortex effects, the linearized forces and moments 
are corrected to reflect losses in suction and the associated formation of 
vortex forces for leading and side edges. The corrections are applied to the 
standard lift, side force and drag coefficients. The corresponding increments 
in the total moment coefficients are calculated by applying the above force 
increments at the appropriate X,Y,Z coordinates for the leading edge stations 
and center of pressure for the side edges. 

For leading edge force calciilations, the lost suction force for each span 
station is given by 


Cg c As'(l -Ks) 


idiere is the coefficient of leading edge suction, C is the local chord, AS'' 
is the local span station width and Ksis the leading edge suction recovery factor. 
This force is subtracted from the direction normal to the section leading edge and 
re-entered as a force component rotated ± 90 ° about the leading edge. The sign 
of the rotation is determined by the sign of the coefficient Ao in the equation 
for leading edge suction. 

The change in the total lift, side force and drag is calculated for each span 
station and is written as a function of four coordinate system rotations whose 
rotation angles are loiown from the leading edge geometry. The origin of each 
coordinate system is located at the leading edge of the section camber line. 

The first transformation involves the rotation of the system (X4,Y4,Z4), whose 
X axis is tangent to the local normal camber line, to the system (X3, Y3, Z3) , 
whose X axis is tangent to the corresponding chord plane as indicated in the following 
sketch: 



94 



where ^ r: JA {[ (clX/dx)^+ (dLz/oLx\ +■ (dT/cLx) ']/'C0SJL'J 

''F 


and 


Cd.-z./d.x')^ 

(cCZ/dLx) 

®F 

-A. 


is streamwise slope due to camber 
is streamwise slope due to twist 
is streamwise slope due to flap deflection 

is the local leading edge sweep angle. 


The sweep terra converts the total streamwise slope to a slope measured in the 
direction normal to the leading edge. 

The two coordinates systems are related by the following transformation 
matrix: 



The second transformation involves the rotation of the system CX3, Y3, Z3), 
whose Y axis is tangent to the leading edge, to the system (X2, Y2, Z2) , whose Y 
axis is normal to the configuration center line and in the plane of the surface. 
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The two coordinates systems are related by the following transformation 
matrix: 









COS JL 


o 

Cy 

> = 

“SlNJL 

COSJL 

o 



0 

o 

1 



•I* 


_ 







The third transformation involves the rotation of the system (X2, Y2, Z2), whose 
Z axis is normal to the local surface plane, to the sytem (Xi, Yi, Zi), vihose X, Y^t 
and Z axes are in the body axes direction. 



The rotation is about the X£, X^ axis and of magnitude O, the local dihedral angle. 
The two coordinate systems are related by the following transformation matrix: 



o o 


• 

Cx 

COSG SIHQ 


Cy 

'z 

SIH© COS© 







The fourth and final transformation involves the rotation of the body axis 
system (Xx, Yx, Zx) to the wind axes system CD, Y, L) . 
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The rotation is about the Y^, Y axis and of magnitude oc , the angle of attack. 
The two coordinate systems are related by the following trans format ion matrix: 


«. .. 

r 


COS ^ 

0 

SlN«k 



0 

-< 


0 

1 

0 

< 

>- 

0 




0 

QOS ok 




The conposite transformation between the (X4, Y4, Z4) coordinate system and the 
(D, Y, L) coordinate system can then be expressed as 


Cd 

■< Cy 


_a 


- 

Cx 
Cy \ 

C-u 






where Sh is the rotation matrix obtained from iraaltiplication of the four previously 
specified transformation matrices. 

Expressing CX4, CY4 and CZ4, in terms of the leading edge suction parametors, 

= C .5 c 

Cy = O 

As (I - Ks) 

we can now write the change in drag, side force and lift resulting from the force 
rotation at each span station: 

= C,c As'O -Ks) Jlp 

A c^ = CjCAs'C' 

A = c., C- As'(l-Ks)il_^ 
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where 


= ' j^cos oC (cos^cos S) -t SIN^(-SIK© SlW-A. Cos S + CosesiHS)] 

4- [cosoc (-COSJL SIN 5") V SJN-c^( SIN© SlHJL SlHS -h Cos © cos 


_r> = [cos© sinJLCoss: 4- 

+ ^ f-C-OS © SIN-A. ^ ■hSlW©COSS'l 


and 


— [-SlNo<*(^OS -A. COS + COSo^(-SJW©SINJCCoS 6 + COS©SlM 

4- (-CosJL Sll^ r) + COS<<( SIN© SIM JA SIN S' 4-Cos ©CoS^)J 

1A.| 


For side edge force calculations, the lost suction force at each chord station 
is give by 

AJL 
^ c 

where Cg is the coefficient of side edge suction, Cy is the tip chord and is 

the local nondimens ional chord increment over which Cs is acting. This force is 
siabtracted from the direction normal to the tip chord and re-entered as a force 
coiqjonent rotated ±90° about the tip chord. Ihe sign of the rotation is 
deteimined by the sign of the coefficient Coq in the equation for side edge suction. 

In a manner similar to that for the leading edge forces, the change in the total 
lift, side force and drag coefficients is calculated for each chord increment and is 
written as a function of three coordinate system rotations whose angles are known 
from the tip geometry. The origin of each coordinate system is located on the chord 
line at the beginning of each chord increment. 
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The first transformation involves the rotation of the system (X3, Y3, Z3) 
whose X axis is parallel to the local camber line, to the system (X2, Y2, Z2) , 
whose X axis is tangent to the tip chord. 



vfliere r= TA N"' f (d-L/d-xX + Ut/'d.xX + (dx/dx) 3 

'• « Op * 


and 


(d.x/ax)^ 
(d.z/ol y)g 


is streamwise slope due to camber 
is streamwise slope due to twist 


(az/ix)y^ 


is streamwise slope due to flap deflection 


The two coordinate systems are related by the following transformation matrix: 




COSS 

0 

-smy 


'' 

^x. 



0 

I 

0 





SiM? 

0 

cos S' 




The second transformation involves the rotation of the system (X2, Y2, Z2), 
vdiose Y axis is normal to the tip chord, to the system (Xi, Yi, Zi), whose 
X, Y, and Z axes are in the body axes direction. 
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The rotation is about the X2 , Xi axis and of magnitude & , the local dihedral 
angle. The two coordinate systems are related by the following transformation 
matrix: 


^x, 





I O O 

o cos & ~SlUO 

o SI we cos© 




Cx. 


The third and final transformation involves the rotation of the body axes 
system (X^, Zi) to the mnd axes system (D, Y, L) . 

L 


D 

X, 

The rotation is about the Yi, Y axis and of magitude , the angle of attack. 
The two coordinate systems are related by the following transformation matrix: 






COSo^ 0 




-< 

Cv 


0 1 

0 




0 

r 


0 

CJOSo<- 


I""’* 


L • J 


— 

— 


L J 

The transformation between the (X3, Y3, 

Z3} coordinate 

coordinate system 

can then be expressed as 





r 


— 



^ - 







Cx, 

-< 

0 

-< 

>- 

f "^“1 














y 


where T is the rotation matrix obtained from multiplication of the three previously 
specified transformation matrices. 
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Expressing CX3, Cyj and C23 

in terras of the side edge suction parametors. 

II 

u 

0 

^3 


- 

- C ^ A K 

= 

Sr:, c A X. 

2 

K1 ^ 

we can now write the change in drag, side force and lift resulting from the force 
rotation at each side edge station: 


= c c ' AX. T 


T c ' D 



AC, 


where 

= T SiH e 

• D 

c«. r-cos**. SINS +SINo^C 0 S© cos s 


T. = 
T, = 


T cos s + SIN a COS S’ 

TCOS<SINe+ 5)22. fsiNol. SIN ? + CoSot C0S9 cos S 1 


The minus sign on the first term of each equation is for the right side of the 
configuration and the positive sign is for the left side. These force increments 
are numerically integrated along each tip chord to obtain the total change in 
lift, side force and drag due to side edge force roation. 
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APPENDIX C 


HYPERSONIC FINITE ELEMENT ANALYSIS 


High ftfach niimber analysis has a number o£ optional methods for calculating the 
pressure coefficient. In each method the only geometric parameter required is 
the element inpact angle, 6, or the change in the angle of an element from a 
previous point. 

The methods to be used in calculating the pressure an impact (6 > 0) and shadow 
(6 < 0) regions may be specified independently. A summary of the program pressure 
options is presented below. 


Impact Flow 


1. Modified Newtonian 1. 

2. Modified Newtonian+Prandtl- Meyer 2. 

3. Tangent wedge 3. 

4. Tangent-wedge empirical 4. 

5. Tangent-cone empirical 5. 

6. OSU blunt body empirical 6. 

7. Van Dyke Unified 7. 

8. Blxont-body shear force 8. 

9. Shock-expansion 9. 


10. Free molecular flow 

11. Input pressure coefficient 

12. Hankey flat -surface empirical 

13. Delta wing empirical 

14. Dahlem-Buck empirical 

15. Blast wave 

16. Modified tangent- cone 


Shadow Flow 

Newtonian (Cp = 0) 

Modified Newtonian+Prandtl -Meyer 

Prandtl-Meyer from free-stream 

OSU blunt body empirical 

Van Dyke Unified 

High Mach base pressure 

Sh o ck - e xp an s i on 

Input pressure coefficient 

Free molecular flow 


A brief review of these methods will be presented in the following text.. 


MODIFIED NEWTONIAN 

This method is probably the most widely used of all the hypersonic 
force analysis techniques. The major reason for this is its simplicity. 
Like all the force calculation methods, however, its validity in any 
particular application depends upon the flight condition and the shape 
of the vehicle or component being considered. Its most general ap- 
plication is for blunt shapes at high hypersonic speed. The usual 
form of the modified Newtonian pressure coefficient is 

2 

C = K sin 6 
P 
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In true Newtonian flow (M =®, V= 1) the parameter K is taken as 2. 
In the various forms of modified Newtonian theory, K is given values 
other than 2 depending on the type of modified Newtonian theory used. 
K is frequently taken as being equal to the stagnation pressure co- 
efficient. In other forms it is determined by the following relation- 
ship (Reference 19 ). 


K 



where 



nose 


6 

nose 


the exact value of the pressure 
coefficient at the nose or leading 
edge 

impact angle at the nose or leading 
edge 


In other work K is determined purely on an empirical basis. 


K = fn (M, a, shape) 

When modified Nev/tonian theory is used, the pressure coefficient in 
shadow regions (6 is negative) is usually set equal to zero. 

MODIFIED NEWTONIAN PLUS PRANDTL-MEYER 

This method, described as the blunt body Newtonian + Prandtl- Meyer 
technique, is based on the analysis presented by Kaufman in Reference 
20 . The flow model used in this method assumes a blunt body with a 
detached shock, followed by an expansion around the body to supersonic 
conditions. This m.ethod uses a combination of modified Newtonian and 
Prandtl- Meyer expansion theory. Modified Newtonian theory is used 
along the body until a point is reached where both the pressure and the 
pressure gradients match those that would be calculated by a continuing 
Prandtl- Meyer expansion. 

The calculation procedure derived for determining the pressure co- 
efficient using the blunt body Newtonian + Prandtl-Meyer technique 
is outlined below. 


1. Calculate free-stream static to stagnation pressure ratio 

V 

r 1y- 1 

2 

p - ffl r 

Po (Y + 1) M^^ 

2. Assume a starting value of the matching Mach number, M 
(for Y= 1.4 assume = 1.35) 
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3. 


Calculate matching point to free-stream static pressure ratio 


Q 



2 

2 + (Y - 1) M ^ 



4. Calculate new free-stream static to stagnation pressure ratio 


P 


c 


Q 1 


M Q 

2 

- 1) (1 - Q) 


5, Assume a new matching point Mach number (1.75) and repeat 
the above steps to obtain a second set of data. 

6. With the above two tries use a linear interpolation equation to 
estimate a new matching point Mach number. This process is 
repeated until the solution converges. 


7. 


Calculate the surface slope at the matching point 


sin 6 

q 


Q - P 
1 - P 


8. Use the Prandtl- Meyer expansion equations to find the Mach 

number on the surface element, M ^ 

o 

• 

9. Calculate the surface pressure ratio 




Y 

Y- 1 


where 

is provided as an empirical correction factor 

Pg is the pressure on the elem.ent of interest 

10. Calculate the surface to free-stream pressure ratio 
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11 . 


Calculate the surface pressure coefficient 




The results of typical caJculations using the above procedure are 
shown in Figure 1. Note that the calculations give a positive pres- 
sure coefficient at a zero impact angle. As pointed out in several 
references these results correlate well with test data for blunt 
shapes. However, if the surface curvature changes gradually to 
zero slope some distance from the blunt stagnation point the pres- 
sure calculated by this method will be too high. This is caused by 
characteristics near the nose intersecting the curved shock system 
and being reflected back onto the body. If the zero slope is reached 
near the nose (such as in a hemisphere or a cylinder) this effect has 
not had time to occur. 

TANGENT-WEDGE 

The tangent-wedge and tangent-cone theories are frequently used to 
calculate the pressures on two-dimensional bodies and bodies of 
revolution, respectively. These methods are really empirical in 
nature since they have no firm theoretical basis. They are suggested, 
however, by the results of more exact theories that show that the 
pressure on a surface in impact flow is primarily a function of the 
local impact angle. In this program the tangent-wedge pressures are 
calculated using the oblique shock relationships of NACA TR-1135 
(Reference 21), The basic equation used is the cubic given by 

(sin^ 0 + b(sin^e + c (sin^ 0 ) +d = 0 or 

' s' s' s 

R^ + bR^+cR + d =0’ 

where 

0 

s 

6 


b 


c 


d 


shock angle 
wedge angle 


M + 2 V • 2 , 

— - V sin o 


M 


2 M + 1 


+ 


(V + 1) . V - 1 

A ' •"» 


• 2 ; 
sin 0 


M 


2. 

cos o 


M' 
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Figure 1, Blunt Body Newtonian + Prandtl-Meyer 
Pressure Results 
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The roots of the above cubic equation may be obtained by using the 
trigonometric solution procedure (see Reference 22) as indicated 
below. 


II 

2 -J- g/3 

cos (w/3) - b/3 


^2 = 

-2 y- p/3 

cos (w/3 + 60°) 

- b/3 

73 = 

-2 V- p/3 

cos (w/3 - 60°) 

- b/3 


II 

7l - b/3 



^^2 

72 - b/3 



R 3 

73 - b/3 



where 

y^ = roots of the reduced cubic equation 

b" , 

p = _ _ + c 

q = 2(b/3)^ - ^ + d 

q 

cos u) = - --.n..— ^ 

2 V-(p/3)-^ 

2 

Rj = sin 6 ^ = roots of the cubic equation 

The smallest of the three roots corresponds to a decrease in entropy 
and is disregarded. The largest root is also disregarded since it 
never appears in physical actuality. 
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For small deflections, the cubic solution becomes very sensitive to 
numerical accuracy; that is, to the number of significant digits car- 
ried. Since this is dependent on the pa.rticular machine emplo , ed, 
an alternate procedure is used. 


When the flow deflection angle is equal to or less than 2,0 degrees, 
the following equation is used instead of the above cubic relationships 
(Reference 23)': 


2 

sin 0 
s 



, y+ 1 d 


Once the shock angle is obtained the remaining flow properties may 
be found from the relationships of Reference 21. 


density = ?2 


6 sin^ 0 
s 

sin^ 6 +5 
s 


temperature = T 2 


7(M^ sin^ e - 1) (M^ sin^ 0, + 5) 

S >5 

36 sin^ 0 

s 


pressure 

coefficient 


S 


7M^ sin^ 0-1 
s 

6 

0.7 


where 


( )2 = conditions behind the shock 


Oblique shock detachment conditions are reached when no solution 
may be found to the above cubic relationships. Under these conditions 
the program uses the Newtonian + Prandtl- Meyer method for continued 
calculations. 
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TANGENT-IVEDGE, TANGENT -CONE, AND DELTA WING 
neWtonian'empirical method 


The tangent -cone and the tangent-wedge Newtonian empirical methods 
used in this program are based on the empirical relationships derived 
below. 





For wedge flow 


sin 6 


sin 0g= ^ 


w 


€) cos (6 - 6 ) 

' ' s w' 


where 


y - 1 

y + 1 


1 + 


(Y - 1) M 


ns 


For cone flow (thin shock layer assumption) 


sin 5 


sin0 = 
s 


'■ (1 - T ) cos (0 - 6 ) 

2 ' ' s c' 


Y - 1 

In the limit as M— n, ^ and cos (0 - 6) = 1 

lim Y+ 1 ' s 


Therefore 


sin0 = 


wedge 
Y + 1 


sin 6 sin 0 = 


w 


cone 


2(Y + 1) . ^ 

Y + 3 
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These limiting expressions for 0 may now be compared with the 
data of TR-1135 (Reference 21) at Y = 7/5 using the following 
similarity parameters. The exact equations contain three vari- 
ables — 0s> 6, and € . Noting that for Y = constant, e = fn(Mj^g) 
only, the preceding equations may be rewritten in the following 
form: 


wedge 


¥ 


cone 


M sin 6 

» j _ w 

^^ns (1 - e) cos (6 - & ) 

s w' 


M 


ns 


M sin 6 

c 

(1 - |) cos (e^ - 6^) 


The parameter (0 - 6) is approximately constant and independent of 
M except near the shock detachment condition. The equations es- 
sentially contain only two variables, M and M sin 6. These are 
used as coordinates to plot the data for wedge flow shown in Figure 
2, A similar plot could be obtained for cone flow. From the figure 
it is seen that the data are nearly normalized with the use of these 
coordinates. 


For rapid calculations we need relationships for as a function of 

M sin 6 that satisfy the following requirements: 

1. The effect of shock detachment is neglected 

2. At M sin 5 = 0, M =1 

ns 

3. The solution asymptotically approaches the M = co line 

d M 

4. Have the correct slope, — : — ?at M sin 6 = 0 

^ ’ d M sin 5 


These conditions lead to equations of the following form 

K 


w 


M' 


wedge M 


ns 


cone M 


ns 
where 


K M' + e 
w 


K 


Y + 1 


w 


K M‘ + e 
c 


■K M* 


M' = M sin 6 

= 2 (Y + 1)/(Y + 3) 
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These expressions are compared with the data of TR-1135 in Figures 
3 and 4 . The cone data are also shown in Figure 5 with the same 
scales as in Figure 2 . 

The pressure coefficient may now be obtained by the following 
relationships for a wedge and cone respectively. 



-1) / m2 


2 sin^6 


1 


(Y- 1) + 2 

4(Y+1) M^32 


Experimental resTilts have shown the pressure on the centerline of a delta 
wing to be in agreement with two-dimensional theory at small values of 
the similarity parameter (M' < 3.0) and with conical flow theory at higher 
values. The previous expressions derived for wedge and cone flows have 
been combined to give these features. The resulting relationships are 
given below. 

K 

-{K^ - -^) M' 

M = K_ M‘ + e ^ 

ns C 


For V = 7/5 


M 


ns 


1 . 09 M sin 5 + e 


-0. 49 M sin 6 


The similarity parameter relationship for pressure is 

“"S = (Yfl),(“n/ - ') 

The shock angle and pressure coefficient calculated from the above 
equations are compared with the experimental results (Reference 28) 
in Figures 6 and 7 , respectively. 
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e 3. Wedge Flow Shock Angle Empirical Correlation 
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Figure 4. Conical Flow Shock Angle Enapirical Correlation 
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Figure 5 . Conical Flow Shock Angle Empirical Correlation 
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Figure 6- Delta W-ing Centerline Shock Angle Correlation 
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Figure 7. Delta Wing Centerline Pressure Coefficient Correlation 
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OSU BLUNT BODY EMPIRICAL METHOD 


The OSU (Ohio State University) blunt body empirical equation 
describes the pressure distribution about cylinders in supersonic 
flow. The equation was presented in Reference 25) and was stated 
to match "all the data obtained on the cylinders in the present test 
series with a maximum deviation of 2. 5 percent. " The expression used is 



0.32 + 0.455 cos6 + 0.195 cos 20 + 0.035 cos 36 - 0.005 cos 46 


where 


0 = peripheral angle on a cylinder 

(= 0 at the stagnation point) = 90® - 6 

Pj = surface pressure 

p = total pressure rise through normal shock 

^ CO 


The pressure coefficient is calculated from the relationship 


C 

P 


where 




) 


t y , 

— = + 1 
P 2 

^co 


K 


Y 


stagnation pressure coefficient = C 

^stag 


freestream pressure 
ratio of specific heats =1.4 


VAN DYKE UNIFIED METHOD 

This force calculation method is based on the unified supersonic- 
hypersonic small disturbance theory proposed by Van Dyke in 
Reference 26 as applied to basic hypersonic similarity results. 

The method is useful for thin profile shapes and as the name implies 
extends down to the supersonic speed region. 
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The similarity equations that form the basis of this method are 
derived by manipulating the oblique shock relations for hypersonic 
flow. The basic derivations are shown on pages 753 and 754 of 
Reference 31. The result obtained for a compression surface under 
the assumption of a small deflection angle and large Mach number 
is (hypersonic similarity equation). 



where H is the hypersonic similarity parameter given by M5. The 
contribution by Van Dyke in Reference 26 suggests that this relation- 
ship will also be valid in the realm of supersonic linear theory if 
the hypersonic similarity paramete r M6 is replaced by the tinified 

supersonic-hypersonic parameter ('/ m^ - l)o. This latter param- 
eter is used in the calculations for this force option in the arbitrary 
body program. 


A similar method may also be obtained for a surface in expansion 
flow with no leading edge shock such as on the upper side of an air- 
foil. The resulting equation is 





V- 1 


H) 


2Y 

Y - 1 



where again H 


is taken to be 


(\/" VI^ - l) 


6 in the 


unified theory approach. 


SHOCK-EXPANSION MEIHOD 

This force calculation method is based on classical shock-expansion 
theory (see Reference 27 In this method the surface elements are 
handled in a, "strip-theory” manner. The characteristics of the first 
element of each longitudinal strip of elements may be calculated by 
oblique shock theory, by conical flow theory, or by a Prandtl-Meyer 
expansion. Downstream of this initial element the forces are cal- 
culated by a Prandtl-Meyer expansion. 

By a proper selection of the element orientation the msthod may be 
used for both wing-like shapes and for more complex body shapes. 

In this latter case the method operates in a hypersonic shock-expansion 
theory miode. 


119 



FREE MDLECULAR FLOW METHOD 


At very high altitudes conventional continuum flow theories fail and 
one must begin to consider the general macroscopic mass, force, 
and energy transfer problem at the body surface. This condition 
occurs when the air is sufficiently rarefied so that the mean free 
path of the molecules is much greater than a characteristic body 
dimension. This condition is known as free molecular flow and the 
method of analysis selected for this program is described in 
Reference 28. This method was also used in Reference 29. The 
equations used were taken from these references and are presented 
below. 

Pressure coefficient 


C 

P 



2 -f 

■ — S sin 6 

VTT 




-fS sin 5) 
e ' 


2 


+ 


(2 -y(s^ 



Shear force coefficient 


C 


f 


(cos 6)f^ 
\/ttS 


^-(Ssin6) 


1 + erf(S sin 6) 


where 

S 

f 

n 


5 



erf 


speed ratio = \/v/2 

normal momentum accommodation coefficient (-1. 0 for 
Newtonian and = 0. for completely diffuse reflection) 

impact angle 

body temperature, °K 

free -stream temperature, “K 


error function erf (x) = 




dx 


f = tahgential momentum accommodation coefficient 

(=0. for Newtonian flow and 1.0 for completely 
diffuse reflection) 
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The pressure force acts perpendicular to the surface and this direction is 
readily obtained since the element normal has already been determined m 
the geometry subroutines. The shear force acts in the direction of the 
tangential velocity component on the surface and this direction is deter- 
mined by taking successive vector products as follows. 


The procedure is illustrated in the 
accompaning sketch where the in- 
cident velocity vector is defined as 


V = Vx i + Vy j +• Vx k 


and t^e surface normal as 


N = Nx i + Ny j N2 k 



First, a su.rface tangent vector (T) is defined by the cross product of the 
normal and velocity vectors; 


T = Tj^i+Tyj + T^l^ 
where 

Tx = Nx - N^ Vx 


Ty = Vx - Nx 


T = N V - N V 
Z Y X 


Then the direction of the shear force (S) is given by the cross product of 
the surface tangent and normal vectors; 


S = Sx i + Sy j + k 
where 


s = 

T„ 

N., - 

T., 

N,, 

X 

Y 

Z 

Z 

Y 

s„ = 

T., 

N„ - 

T,. 

N.^ 

Y 

Z 

X 

X 

Z 


Tx 

Ny - 

T 

-V 

X 

TVJ 

*^x 
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The final components of the shear force in the vehicle axis system are 
given by 

SHEAR^ 

SHEARy 

SHEAR 2 

where 

SHEAR is the shear force as calculated by the free molecular flow 
equations . 

STOTAL = ^S^^ + Sy^ + S^^ 

In using the free molecular flow method the above analysis must be 
carried out over the entire surface of the shape including the base, 
shadow regions, etc. When the free molecular flow method is 
selected, it is used for both impact and shadow region. 


= (SHEAR) (s^) / STOTAL 
= (SHEAR) (Sy) / STOTAL 
= (SHEAR) (S^) / STOTAL 


The plane foiroed by the velocity vector and the surface normal is 
referred to as the velocity plane (shaded region in the sketch), 
since both the incident and surface velocity are in this plane. This 
definition is correct for two-dimensional flow, however, it is only 
an approximation to the shear direction in the general arbitrary -body- 
case . 


HANKEY FLAT-SURFACE EMPIRICAL METHOD 

This method uses an empirical correlation for lower surface pressures on 
blunted flat plates. The method, derived in Reference 30, approximates 
tangent-wedge at low impact angles and approaches Newtonian at high impact 
angles. The pressure coefficient is given by 



nMIT.HM-BIICK EMPIRICAL METtOJ) 


This is an impact method that has been derived such that tangent-cone and 
Newtonian results are approximated, respectively, at low and high values 
of the impact angle. The empirical relationships presented in Reference 31 


are 

for 

6 < 22. 5° 

O 

■P 

II 

1 + (sin , . „,5/4 

3/4 (sin 6) 

(4 cos 6 cos 26) 

for 

6 > 22.5° 


-,0-2 
2. 0 sin 5 


BLAST WAVE PRESSURE INCREMENTS 


This method uses conventional blast-wave parameters to calculate the over- 
pressure due to bluntness effects. Force contributions determined by this 
procedure must be added to the regular inviscid pressure forces (tangent- 
wedge, tangent-cone, Newtonian, etc. ) calculated over the same vehicle 
geometry. The specific blast wave solutions used in the Program were de- 
rived by Lukasiewicz in Reference 32: 


where 


= A M„ 


(Cd) 


Rt j... 


(Xq - X)/d 


+ B 


Cjj is the nose drag coefficient 
d is the nose diameter or thickness 
Xq is a coordinate reference point 


and the coefficients A, B are 


Flow 

j 

A 

B 

Two-Dimensional 

0 

0.121 

0.56 

,Axisymmetric 

1 

0.067 

0.44 
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MODIFIED TANGENT-CONE METHOD 

This method, originally developed for use on cones with elliptical cross- 
sections, modifies the tangent-cone result by an increment representing 
the deviation from an average pressure divided by an average Mach number. 
More specifically, the following equations are used (after Jacobs, Reference 
33): 


■ ^Ptc 


^Ptc ~ ^Pavg 


M 


avg 


where 


C is the surface pressure coefficient 
P 

is the conventional tangent-cone pressure coefficient 
Cpg^^g is the average pressure coefficient 
= ZCp^ a/£a , A is element area 


M 


avg 


is the average Mach number, defined for an equivalent 
cone having pressure coefficient Cp^^^. 


HIGH MACH BASE PRESSURES 

For a body in high speed flow it might be expected that any base regions 
would experience total vacuum. That is. 



1 


— 

2 ® 


However, the viscosity of real gases causes some pressure to be felt in 
base region and experimental data have shown this to be roughly 70% 
vacuum for air. Therefore, the expression 



has been included in the program. 
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